The Yield-Strain in Shear Banding Amorphous Solids 
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In recent research it was found that the fundamental shear-localizing instability of amorphous 
solids under external strain, which eventually results in a shear band and failure, consists of a 
highly correlated array of Eshelby quadrupoles all having the same orientation and some density 
p. In this paper we calculate analytically the energy E(p, 7) associated with such highly correlated 
structures as a function of the density p and the external strain 7. We show that when the strain 7 is 
smaller than some characteristic 7 Y the minimum energy solution is attained for p — (i.e. isolated 
localized plastic events). For 7 > 7 Y there is a bifurcation allowing a finite density of quadrupoles. 
The paper culminates with an analytic estimate of 7 Y , without any free parameters, in terms of 
material characteristic constants. 



I. INTRODUCTION 

Amorphous solids are obtained when a glass-former is 
cooled below the glass transition [l]-[3| to a state which 
on the one hand is amorphous, exhibiting liquid like or- 
ganization of the constituents (atoms, molecules or poly- 
mers), and on the other hand is a solid, reacting elasti- 
cally (reversibly) to small strains. There is a large vari- 
ety of experimental examples of such glassy systems, and 
theoretically there are many well studied models 
based on point particles with a variety of inter-particle 
potentials that exhibit stable supercooled liquids phases 
which then solidify to an amorphous solid when cooled 
below the glass transition. Typically all these materials, 
both in the lab and on the computer, exhibit a so-called 
yield-stress above which the material fails to a plastic 
flow. In previous research [7l-fIoj it was pointed out that 
depending on the protocol of cooling to the glass state, 
the plastic response of the system can be either via ho- 
mogeneous flow or via shear bands. The former obtains 
typically when the quenching to the glass state is "fast" , 
whereas the latter when the quench is "slow" . The plastic 
response via a shear band is exemplified in FigJT] which 
exhibits a typical failure of a sample of metallic glass 
when subjected to compressive stress. When the stress 
exceeds some yield stress, the sample, rather than flowing 
homogeneously in a plastic flow, localizes all the shear in 
a plane that is at 45 degrees to the compressive stress 



FIG. 1: A typical example of the failure of a metallic glass 
sample when subjected to compressive stress. The material 
localizes the stress in a plane that is at 45 degrees to the 
compressive stress axis, and then breaks along this plane 



axis, and then breaks along this plane [l[. 



In recent work |T2j it was argued that the fundamental 
instability that gives rise to shear bands is the appearance 
of highly correlated lines of Eshelby quadrupoles (and 
see below for precise definition) which organize the non- 
affine displacement field of an amorphous solid such that 
the shear is highly localized along a narrow band. How 
this fundamental instability results in macroscopic shear 
bands, why these appear in 45 degrees to the principal 
stress axis, and what determines the difference in plastic 
response between fastly and slowly quenched glasses are 
all subjects of this paper. We will also present an ab- 
initio calculation of the yield stress at which an amor- 
phous solid is expected to response plastically with shear 
localization. 



In Sect. [TT] we review briefly the type of numerical 
simulations that we do and explain the basic facts about 
plasticity of amorphous solids. Section IIIII exhibits the 
fundamental solution of an Eshelby quadrupolar plastic 
event. This section is not particularly new but is im- 
portant for our purposes in setting the stage and the 
notation for the next Sect. IIVI in which we compute 
the energy of N such Eshelby quadrupoles in the elastic 
medium. We show explicitly that as a function of the ex- 
ternal strain (or the resulting stress) there is a threshold 
value at which a bifurcation occurs. Below this value only 
isolated Eshelby quadrupoles can appear in the system, 
leading to localized plastic events. Above this thresh- 
old a density of such quadrupoles can appear, and when 
they do appear they are highly correlated, preferring to 
organize on a line at 45 degrees to the principal shear 
direction. In Sect. |V] we present the analytic estimate 
of the yield strain, and demonstrate a satisfactory agree- 
ment with the numerical simulations. Finally, in Sect. 
IVII we provide a summary of the most important results 
of the paper and offer a discussion of the road ahead. 
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II. PLASTICITY IN AMORPHOUS SOLIDS 
AND SIMULATIONS 

As a background to the calculations in this paper we 
need to briefly review recent progress in understanding 
plasticity in amorphous solids [jj H, Il3l - [l5| . Below we 
deal with 2-dimensional systems composed of N point 
particles in an area A, characterized by a total energy 
U (ri, T2, ■ • • r n ) where rj is the position of the z'th parti- 
cle. Generalization to 3-dimensional systems is straight- 
forward if somewhat technical. The fundamental plastic 
instability is most cleanly described in athermal (T = 0) 
and quasi-static (AQS) conditions when an amorphous 
solid is subjected to quasi-static strain, allowing the sys- 
tem to regain mechanical equilibrium after every differen- 
tial stress increase. Higher temperatures and finite strain 
rates introduce fluctuations and lack of mechanical equi- 
librium which cloud the fundamental physics of plastic 
instabilities with unnecessary details [131 ] . 

In our AQS numerical simulations we use a 50 — 50 
binary Lcnard- Jones mixture to simulate the shear local- 
ization discussed in this work. The potential energy for 
a pair of particles labelled i and j has the form 
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where the parameters Aq, A\ and Ai are added to smooth 
the potential at a scaled cut-off of r/a = 2.5 (up to the 
second derivative). The parameters oaa, &bb and gab 
were chosen to be 2sin(7r/10), 2sin(7r/5) and 1 respec- 
tively and 6aa = £bb — 0.5,6,4b = l(see [l6[). The 
particle masses were taken to be equal. The samples 
were prepared using high-temperature equilibration fol- 
lowed by a quench to zero temperature (T = 0.001) (see 
[ijj )■ For shearing, the usual athermal-quasistatic shear 
protocol was followed where each step comprises of an 
affine shift followed by an non-affine displacement using 
conjugate gradient minimization. The simulations were 
conducted in two dimensions (2d) and employed Lees- 
Edwards periodic boundary conditions. Samples were 
generated with quench rates ranging from 3.2 x 10~ 6 
to 3.2 x 10~ 2 ( in LJ units ), and were strained to 
greater than 100 percent. Simulations were performed on 
system-sizes ranging from 5000 to 20000 particles with a 
fixed density of p = 0.976 (in LJ units). The simula- 
tions reported in the paper have 10000 particles and a 
quench-rate of 6.4 x 10 -6 (in LJ units). 

We choose to develop the theory for the case of exter- 
nal simple shear since then the strain tensor is traceless, 
simplifying some of the theoretical expressions. Apply- 
ing an external shear, one discovers that the response of 
an amorphous solids to a small increase in the external 
shear strain Sj (we drop tensorial indices for simplicity) 
is composed of two contributions. The first is the affine 
response which simply follows the imposed shear, such 




FIG. 2: (Color Online). Left panel: the localization of the 
non-affine displacement onto a quadrupolar structure which 
is modeled by an Eshelby inclusion, see right panel. Right 
panel: the displacement field associated with a single Eshelby 
circular inclusion of radius a, see text. The best fit parameters 
are a « 2.5 and e* ~ 0.1. To remove the effect of boundary 
conditions, the best fit is generated on a smaller box of size 
{x, y) € [25.30, 75.92] 



that the particles positions = Xi , yi change via 

$1 Vi = x i 
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This affine response results in nonzero forces between the 
particles (in an amorphous solid) and these are relaxed 
by the non-affine response it; which returns the system 
to mechanical equilibrium. Thus in total — > r[ + u,. 
The nonaffinc response Ui solves an exact (and model 
independent) differential equation of the form [f| [TT| 



d 7 ~ 13 " 3 
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d 2 U( ri ,r 2 , 



— is the so-called Hessian ma- 
is known as the non-affine 



where H^j \ 
trix and = 

force. The inverse of the Hessian matrix is evaluated af- 
ter the removal of any Goldstonc modes (if they exist). A 
plastic event occurs when a nonzero eigenvalue Xp of H 
tends to zero at some strain value jp. It was proven that 
this occurs universally via a saddle node bifurcation such 
that Xp tends to zero like Xp ~ ^/■yp — 7 15]. For values 
of the stress which are below the yield stress the plastic 
instability is seen Q as a localization of the eigenfunction 
of H denoted as which is associated with the eigen- 
value Xp, (see Fig. [2] left panel). While at 7 = all the 
eigenfunctions associated with low-lying eigenvalues are 
delocalizcd, "Jp localizes as 7 — >• 7p (when Xp —> 0) on a 
quadrupolar structure as seen in Fig. [2] left panel for the 
non-affine displacement field when the plastic instability 
is approached. These simple plastic instabilities involve 
the motion of a relatively small number of particles (say 
20 to 30 particles) but the stress field that is released has 
a long tail. 

When the strain increases beyond some yield strain, 
the nature of the plastic instabilities can change in a fun- 
damental way [lOl ] . The main analytic calculation that 
is reported in Sect. IIVI shows that when the stress built 
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FIG. 3: (Color online). The energy vs. strain in a typical nu- 
merical simulation in our system. Shown are the points on the 
curve for a regular plastic event involving a single quadrupo- 
lar structure, (marked in a red circle on the main graph and 
in insert (a) which is a blow up of the same graph), and also 
the point on the curve that results in a plastic instability lead- 
ing to a shear band, yellow ellipse on the main graph and the 
blow up in insert (b). The regular plastic event is not even 
seen without a blow up. 



in the system is sufficiently large, instead of the eigen- 
function localizing on a single quadrupolar structure, it 
can now localize on a series of Af such structures, 
which are organized on a line that is at 45 degrees 
to the principal stress axis, with the quadrupo- 
lar structures having a fixed orientation relative 
to the applied shear. In Fig. [3] we show a typical 
stress vs. strain curve of a glassy material, exhibit- 
ing where one find a localized plastic event and when 
the event is sub-extensive, spanning the system with a 
line of quadrupoles. Fig. 2] shows the non-amne field 
that is identical to the eigenfunction which is associated 
with this instability, clearly demonstrating the series of 
quadrupolar structures that are now organizing the flow 
such as to localize the shear in a narrow strip around 
them. This is the fundamental shear banding instability. 
Note that this instability is reminiscent of some chainlike 
structure seen in liquid crystals, arising from the orien- 
tational elastic energy of the anisotropic host fluid fl8| . 
and ferromagnetic chains of particles in strong magnetic 
fields [lj|. The reader should note that the event shown 
in Fig. 2] will move the particles only a tiny amount, and 
it is the repeated instability where many such event hit at 
the same region which results in the catastrophic event 
shown in Fig. [TJ Nevertheless this is the fundamental 
shear localization instability and in the sequel we will 
have to understand why repeated instabilities hit again 
and again in the same region. We should stress that 
this is not inevitable, for samples that are prepared by 
a fast quench these instabilities appear at random places 
adding up to what seems to be a homogeneous flow. For 
further discussion of this point sec Sect. IVI1 





FIG. 4: (Color Online). Left panel: The nonaffine displace- 
ment field associated with a plastic instability that results in 
a shear band. Right panel: the displacement field associated 
with 7 Eshelby inclusions on a line with equal orientation. 
Note that in the left panel the quadrupoles are not precisely 
on a line as a result of the finite boundary conditions and the 
randomness. In the right panel the series of N = 7 Eshelby 
inclusions, each given by Eq. (|18[) and separated by a distance 
of 13.158, using the best fit parameters of Fig. have been 
superimposed to generate the displacement field shown. 



III. DISPLACEMENT IN 2D FOR A CIRCULAR 
INCLUSION 

As said above, the shear localizing instability appears 
only when the stress exceeds a threshold. To explain 
why, we turn now to analysis. As a first step we model 
the quadrupolar stress field which is associated with the 
simple plastic instability as a circular Eshelby inclusion 



A. Circular Inclusion 

We consider a 2d circular inclusion that has been 
strained into an ellipse using an eigenstrain or a stress- 
free strain e*^ which we take to be traceless ie e* 7 = 
[20| . A general expression for such a traceless tensor can 
be written in terms of a unit vector n n and a scalar e* as 



E a/3 = e * (2 n ahp - 8 a p) ■ 



(4) 



We also assume that a homogenous strain acts glob- 
ally (which in our case also triggers the local transforma- 
tion of the inclusion). This strained ellipsoidal inclusion 
both feels a traction exerted by the surrounding elas- 
tic medium resulting in a constrained strain e c a o in the 
inclusion, and itself exerts a traction at the inclusion- 
clastic medium interface resulting in the originally un- 
strained surroundings developing a constrained strain 
field e^(X). 

A fourth-order Eshelby tensor S a /3js can be defined 
which relates the constrained strain in the inclusion e 
to the eigenstrain e* g viz. 
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t a/3 



7(5 ' 



(5) 



Now for an inclusion of arbitrary shape the constrained 
and displacement field u a inside 



strain eig, stress a a Q 
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the inclusion are in general functions of space. For ellip- 
soidal inclusions, however, it was shown by Eshelby \2l\ - 
l23l | that the Eshelby tensor and the constrained stress 
and strain fields inside the inclusion become indepen- 
dent of space. We work here with a circular inclusion 
which is a special case of an ellipse and hence for such an 
inclusion, the Eshelby tensor is [2l| - |23j 



as e* 5 and are traceless. Note that Eq. (J7J) implies 
that for a traceless eigenstrain, the constrained strain 
inside the inclusion viz. e^g is also traceless and thus we 
obtain from Eq. (TIB"]) , 
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where v is the Poisson's ratio. From Eqs. (|5|) and (J6j) , we 
obtain 



1 + 1/4(1-1/) a/3 l + i/ Q/3 
-Z 7 Z 7 

4 (1 + 1/) (1 - 1/) Q/3 i + i/ a0 



T+T c ^ 

(14) 



fc a/3 



4i/-l 



(1-f) 
3-4// , 

4(1-//) 6(1/3 



4// 



8(1-!/) 



For a traceless eigenstrain. 
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(7) 



The total stress, strain and displacement field inside the 
circular inclusion is then given by 
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Here the super-script J indicates the inclusion and the 
<r*0 denotes the eigenstress which is linearly related to 
the eigenstrain by er*^ = C Q/ 3 7 5e* (5 , and which for an 
isotropic clastic medium simplifies further using 



SpS + SagS/Sy) , (9) 



to: 
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where A and \x are the Lame's parameters. One can cither 
choose the two Lame's coefficients or £ and v as the two 
independent material parameters. The relations between 
them are given by 
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The stress in the inclusion can now be written down 
in terms of independent variables using Eq. (|10j) by 



CT a/3 — CafjyS (e^g - 6*5 + €^g) 
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B. Constrained Fields in the Elastic Medium 

In the surrounding elastic medium the stress, strain and 
displacement fields are all explicit function of space and 
can be written 



^(i) = ^(i) + ^ 
<(x) = <(x) + <>(x) 



(15) 



In order to compute the displacement field u c a {X) in 
the isotropic elastic medium we need to solve the Lamc- 
Navier equation 



£ 



2 (1 + v) (1 - 2v) dX a dX 1 2 (1 + v) 3X & dX t 



= 0(16) 



as there are no body forces present in our calculation. 
The constrained fields in the inclusion will supply the 
boundary conditions for the fields in the elastic medium 
at the inclusion boundary. Also as r — > oo the con- 
strained displacement field will vanish. 

All solutions of equation Eq. ([To]) also obey the higher 
order bi-harmonic equation 



dXgdXgdX^dX 



= v 2 v V„ = o. 



(17) 



Thus our objective is to construct from the radial so- 
lutions of the bi-laplacian equation Eq. (Tl~7|) derivatives 
which also satisfy Eq. (fT6"|) . Note that Eq. (fTT)) is only a 
necessary (but not a sufficient) condition for the solutions 
and Eq. ([16]) still needs to be satisfied. The calculation 
is presented in Appendix A, with the final result 
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Fit to the data 



Armed with this analytic expression we return now to 
our numerics, cf. Fig. [5J and fit the two parameters 



5 



in Eq. (|T8|) to the data of the displacement exhibited 
by a single localized plastic event. The result of this 
procedure is a « 2.5 and e* « 0.1. The quality of this 
fit can be judged from the right panel of Fig. [5] where 
we exhibit the form of Eq. (|18|) with the parameters 
fitted to the displacement field in the left panel. Also 
the value of a appears reasonable since it means that 
about na 2 « 20 particles are involved in the core of the 
relaxation event. On physical grounds this is about the 
right order of magnitude. 

We will keep these parameters fixed in all our calcula- 
tions below. The reader should note that this is an ap- 
proximation when there are multiple quadrupoles in the 
system, since they influence each other and the solution 
leading to Eq. (fT8|) should be repeated in the presence of 
many quadrupoles. We expect however that the changes 
in the parameters should not be large when the density 
of the quadrupoles is small. We will always work in the 
small density limit pa 2 <C 1 where p is the area density 
of quadrupoles N / 1 L 2 . 

IV. THE ENERGY OF JV ESHELBY 
INCLUSIONS EMBEDDED IN A MATRIX 

A. Notation 

Having the form of a single quadrupole, Eq. (1181) . 
we turn now to the calculation of the energy associ- 



ated with Af quadrupoles embedded in an clastic ma- 
trix. Once computed, we will show later that for large 
strains the minimum of this energy is obtained for a line 
of quadrupoles all having the same orientation. From 
now on we use the notation that u a .p = ■ The en- 
ergy of Af Eshclby inclusions embedded in a linear clastic 
medium (or matrix) Af, is given by the expression 



»=1 JV J V l-i = l v o 



where the superscript i indicates the index of the inclu- 
sion and m indicates the matrix. We evaluate Eq. Q19p 
in Appendix B. The result can be expressed in terms 
of four contributions: E m &t which is the contribution of 
the strained matrix, E°° which is the energy of the Af 
qudrupolcs in the external strain, i? cs h which represents 
the self energy of the Af quadrupoles (their cost of cre- 
ation) and lastly E- lnc represents the energy of interaction 
between the inclusions. Explicitly 
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The above expressions are specific to 2D, for a global strain corresponding to simple shear under the linear 



approximation. Thus e°ff y — 5r here, and the traceless eigenstrain takes the form — 2e*nx 1 n y ' 



The form of E lnc is not final, and we bring it to its final form in Appendix C. The final result is 



FIG. 5: Schematic of four Eshelby inclusions embedded in a matrix 'm' 




where f y = 

Our task is now to find the configuration of N Eshelby 
quadrupoles that minimize the total energy. Obviously, 
if the external strain 7 is sufficiently large, we need to 
minimize E°° separately, since it is proportional to 7. 
The minimum of (|21[) is obtained for 

n^= n y = -L. (25) 
I 



Substituting this result in Eq. (f24)) simplifies it consid- 
erably. We find 




We can find the minimum energy very easily. Denote B[x — x ]. The minimum is obtained at x = 1/2, or 
x = (h ■ r) 2 , and minimize the expression A[2x — l] 2 — 
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cos0 = y/l/2. We thus conclude that when the line of 
correlated quadrupole forms under shear, this line is in 
45 degrees to the compressive axis, as is indeed seen in 
experiments. Note that for other external strains which 
are not consistent with a traceless strain tensor (or in 
3-dimensions) this conclusion may change. 

The physical meaning of this analytic result is that 
it is cheaper (in energy) for the material to organize Af 
quadrupolar structures on a line of 45 degrees with the 
compressive stress, all having the same orientation, than 
any other arrangement of these Af quadrupoles, includ- 
ing any random distribution. This explains why such 
a highly correlated distribution appears in the strained 
amorphous solid, and why it can only appear when the 
external strain (or the built-up stress) are high enough. 
This fact, in addition to the observation that such an 
arrangement of Eshclby quadrupoles organizes the dis- 
placement field into a localized shear, explains the origin 
of this fundamental instability. 



For an isotropic matrix, using Eq. (|9]), we have 
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for a symmetric traceless eigenstrain. Using Eq. ([281 
circular inclusions each of radius a in 2D, we obtain 
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V. ESTIMATE OF YIELD-STRESS AND 
NUMBER OF ESHELBY QUADRUPOLES 

In this section we turn to estimate the yield stress and 
the associated density of Eshelby quadrupoles. To this 
aim we need to compute one other energy term that was 
not needed until now, namely E os h which was the same 
for all the configurations of the quadrupoles. 



For all eigenstrains equal (ie. e 1 -*' 1 ' = e*), we have 
£ira 2 Af(e*) 2 
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4(1 -i/ 2 ) 



(31) 



A. Expression for E es h 



B. Ei nc for a line of equidistant quadrupoles with 
the same polarization 



The energy term E cs h was given by Eq. 
be re-written as 



2 -V 
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c 7(5 



which can 



At this point we need to compute the energy term E{ nc 
for the special configurations of Af quadrupoles that are 
cqui-distant and with the same polarization, organized 
(27) in a line. In this case we have 
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Cf. Eq. 



n, 



f , and Ri 



\j - i\R 



(32) 



Starting from Eq. (|26|) we specialize to the present situ- 
ation 
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For (n • f) 2 = 1/2, the above expression reduces to 
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We have 



Af-1 Af 
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where Q{s) is the Ricmann zeta function. Thus we obtain 
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At this point we realize that the distance R between 
the quadrupoles is not determined. We will choose R 
by demanding that the area density of the quadrupoles 
remains intensive in the thermodynamic limit, or p = 
N/L 2 -^Const. We thus choose R = L 2 /Na = 1 /{pa) 
and find 
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Similarly from equation 13 1[ we obtain 
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From equation (|25[) we have 
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Thus the plastic energy density is given by 
E(p,j) _ E°° + E csh + E inc 



(37) 



(38) 



(39) 



L 2 ~ 
£ir(e*) 2 



L 2 



(40) 



A ( 1 - 2- ) {pa 2 ) - B{pa 2 Y + C{pa 2 f 



4(1 -1/2) 

where A = 1, B = 4£(2) and C = 6C(4). 7 Y is defined as 



7v 



2(1 - v) 



(41) 



Eq. (j40|) is plotted, using the numerical values of the 
parameters found in our numerical simulations in Fig. 
IV Bl for various values of 7. For 7 < 7 Y the minimum 
of the expression is attained always at p = 0, and the 
shear localization cannot occur. Only at 7 > j Y a new 
solution opens up to allow a finite density of the Es- 
helby quadrupoles. Therefore 7 Y is by definition the yield 
strain. 
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FIG. 6: (Color Online). The total plastic energy for the cre- 
ation of an array of quadrupoles with (areal) density p for 
three values of 7: 7 = j Y — 0.1, 7 = 7 Y — 0.05, and 7 = 7 Y 



VI. SUMMARY AND CONCLUSIONS 

We have presented a theory of the fundamental insta- 
bility that leads to shear localization and eventually to 
shear bands. One remarkable observation is that the nat- 
ural plastic instability that occurs spontaneously in our 
simulations results in a displacement field that is sur- 
prisingly close to the one made by an Eshelby circular 
inclusion, see Fig. [5J The best fit for the parameter a, of 
the order of 2.5 is in agreement with the intuitive belief 
that shear transformation zones involve 20-30 particles, 
as 7ra 2 would predict. Basing our analysis on this simi- 
larity we could develop an analytic theory of the energy 
needed to create Af such inclusions, whether scattered in 
the system randomly or aligned and organized in highly 
correlated shear localized structure. We discover that 
the latter becomes energetically favorable when 7 exceeds 
7 Y = 9( i_ jA ■ In our system v sa 0.215, and with our best 

0.07 which is right on the 



2(1-!/) • 

fit e* w 0.1 we predict 7 Y 
mark as one can see from Fig. [3] 

While we believe that our calculation of the energy of 
quadrupoles is accurate for densities such that pa 2 <C 



9 



-19900 




FIG. 7: (Color Online). Repeated shear localization insta- 
bilities in the case of slow quench. Upper panel: the energy 
vs. strain and the instabilities that were chosen for display 
in the lower panel. The average non afhne displacement field 
after the instabilities that are marked in the upper panel. The 
point to notice is that the instability falls repeatedly on the 
same band, accumulating to a shear band. 



1, the interaction between the quadrupolcs become much 
more involved for higher densities, and we avoided this 
complication. The consequence is that we cannot pre- 
dict a-priory the critical density of our quadrupoles, and 
we leave this interesting issue for future research. An- 
other issue that warrants further study is whether the 
parameters a and e* are material parameters which are 
determined by the small scale structure of the glass, and 
if so, how to estimate them a-priori. Also, are these pa- 
rameters dependent on the way the system is stressed, 
i.e. via shear or via tensile compression etc.. 

Next we need to discuss the difference in behavior of 
glasses that were quenched relatively quickly and those 
that were quenched relatively slowly. In our simulations 
we find that the latter exhibit the fundamental shear 
localization instability again and again along approxi- 
mately the same line, accumulating displacement that 
develops into a shear band. On the other hand the for- 
mer tends every time to have plastic instabilities at dif- 
ferent places, sometime localized and sometime more cor- 
related, and on the average this then appears like a ho- 
mogenous flow. This difference is stressed in Fig. IVII 
where the case of slow quench is exhibited. It appears 
that fastly quenched systems may have plastic instabili- 
ties almost everywhere, and there is no special preference 



for one line or another. Slowly quenched systems are ini- 
tially harder to shear localize, but once it happens in a 
particular (random) line it is likely to repeat again and 
again along the same line. Making these words quantita- 
tive is again an issue left to future research. In particular 
it is interesting to find what is the precise meaning of fast 
and slow quenches, fast and slow compared to what, and 
how this changes the microscopic local structure. 

Finally the effects of temperature and finite strain rates 
on the present mechanism constitute a separate piece of 
work which of course is of the utmost importance. Like 
the other open subject mentioned above, it will be taken 
up in future research. 
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Appendix A: The displacement field of an Eshelby 
circular inclusion 



1. Solutions of the Lame-Navier equation 

We look for linear combinations of derivatives of the ra- 
dial solutions of equation[T7]which are linear in the eigen- 
strain e*^ and go to zero at large radii. In addition the 
terms must transform as components of a vector field. 
Such a solution can be written down as 



AC 



d In? 



Be 



d 3 In \ 



M dX a dX B dX-, 



Ct 



Q3 ( r 2 m r ~J 

M dX a dX B dX. 



(Al) 



where X a is the a component of the position vector with 
the origin at the center of the Eshelby quadrupole, and 
r = \X\. It can be checked that any other terms are 
either zero, do not go to zero as r — > oo, or do not trans- 
form as components of a vector. We re-write equation 
Eq. CH]) as 



1 



ff 2 u c r 



l-2vj dXadXj dXpdXp 



(A2) 



which is the equation for the constrained displacement 
field in the elastic matrix subject to appropriate bound- 
ary conditions. From Eq. (|A1[) . we obtain 
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2 c a 

~ Ae* 



d 2 



dx p dx p 



lni 



VBe 



d 3 



,lX dX a dX v dX x 



i) 2 



In \ 



+ Ce 



g 3 



We need the following identities 
() 2 



dX B dX 



(lnr) = 



if 



(A4) 



(r 2 hi r) = 4 In r + 4 



dXpdXp 
Thus we obtain from equation IA3[ 

d 3 lnr 



dXpdXp 



4Ce 



'/A 



dX a dX ri dX x 



(A5) 



d 2 u c 7 

dxJrxZ 



d 2 



dXadX^ 
d 



Ae* 



d In? 
'~dX~„ 



Be 



,lX dX a dX v dX x 



dXpdXp 
I 

Similarly, the expression for 



d 2 



dx p dx p 



r 2 In \ 



d 3 hn 



d 3 (r 2 In r) 

" A dx~dx„dx x ^ ^ vX dx y dx n dX; 



cc 



dX a 



Ac 



d 2 lnr 

iidx^dXr, 



Be 



a 4 In 



" A dX^dX^dX^Xx 



Ce 



Q4 ^,2 j n , 



dX 7 dX 7 dX n dX x 



(A3) 



(A6) 



which can be re- written (noting that the second and third 
terms involve the laplacian for which we have identities 
from Eqs. (fA4f as 



d 



8X a dX-y 



Ac 



3 2 lnr „ d 2 (4 lnr + 4) 



(a + 4C) £ ; a 



7 " dX 7 dX v 



Ce 



" A ax r ,aA A 



<9 3 In j 



Plugging expressions (| A5|) and (|A7I) in Eq. 
obtain 



2J), we thus 



dX a dX v dX x 



(A7) 
I 



(A+4C) . 8 3 In r 



4Ce! 



3 lnr 



= 



1-2^ C J7A dX a dX n dX x ~ ^^v*dX a dX v dX x 

We can thus re- write Eq. (|Al|) as 

9 In r „ „ <9 3 In r 



A+4C 
l-2v 



The following identities are now required 



+ Be 



AC 



A 



c * d 3 lnr 

"r/A 0X a dX^dX x 



= 0^C = -- 



d 3 (r 2 In r) 

dX a dX fi dX 1 8 (1 - i/) " P7 dX a dX p dX, 



d In r Aa 



9 3 In r 

dxjyx^trx^ 

d 3 (r 2 In r) 
dXadXpdX^ 



— 2r 2 (A Q (5 / 3 7 + Xp5 a ^ + A 7 (5 Q ^) + SAqA^A 



2r 2 {X a 8p 1 + Xp5 ai + A 7 <5 Q/3 ) - 4A Q A (3 A 7 



(1-^) 



(A8) 



(A9) 



(A10) 
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Using these relations we can re- write Eq. IA9I as 



Xfi 



-2r 2 (XaSfr + Xp5 ai + Xy5 a p) + 8X a Xf}X. 



8(1-!/) 
Af* ^ 



2r 2 (X ct (5 j a 7 + X^<5 a7 + X 7 <5 Q( 3) — 4X a XpX 



.4 



2B 



r A 4 (1 - v) r 2 



(X Q (5^ 7 + XpSay + XySafi) 



A 



8B 



_r e 2(\-v)r A 



e*g„X a XpX^AU) 



Remembering that the eigenstrain is traceless,we find that e^ 7 [X a 8^ + Xp5 ai + X 7 (5 Q( g) = 26*^X^9 using which we 
can simplify Eq. (|A11[) to obtain 



A — ML _ A 

r 2 r 4 2(l-i/)r 2 



r 6 1 2(l-i/)r 4 



X a e *p 7 X /3 X -f 



A l-2i/ _ 4B 
2(1-1/) l 71 " 



8B , A 
7^" 2(l-i/)r 4 



X a € tL p 1 XpX 1 



(A12) 



At r = a (the radius of the circular inclusion), the form of expression Eq. (|A12p must match the form of the constrained 
displacement field of the inclusion which from Eq. (J7) is gjfj^ e*^^. Thus the co-efficient of the second term in 
expression (|A12[) must go to zero at the inclusion boundary, which gives us 



8B , A _ n 



B 



16(l-i/) 



(A13) 



Thus we have 



A l-2i/ _ 4 -a 2 A 
T 1 2(l-i/) 16(l-i/) 



8 -q^A , A 
F 5 ' 16(1-1/) 2(l-i/)r 4 



X a e^, y X / 3X 7 



4r 2 (l-i/) 



2(1-21/) + ^ 



£ q/3^ + 2r 4 (l-i/) 



1 - 2_ 



^(2^-^/3 X)> 



(A14) 



And the value of u„ at r = a should match the value obtained from Eq. (JT)), implying 
3-4z/ 



4(1-;/) 4(l-i/) a 2 
The expression for becomes: 



[2 (1 - 2f ) + 1] => — r = T7T ^7 A = a 2 



4(l-i/) 4(1-!/) a 2 



(A15) 



1 /a 2 



a 4(l-i/) Vf 2 
From Eq. HJ we have e*^ = e* (2n a np — <5 a /3) and thus 



1 



2(1-!/) V r 



1-1^ 



; XaXpXy 

E /37 ^2 



(A16) 



2n„ n • X - X 



X P e }-y X l = e * X P ( 2 «/3«7 ~ <W ^7 = 



( fl .x) ! 



(A17) 



allowing us to write the final vectorial expression for the displacement field: 



u c X = 



4(l-i/) V^ 2 



2(1 - 2v) 



2(l-i/) V^ - 



2 nX 



2n(n-X) -X 

2 

X 



1 



(A18) 
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We can also derive expressions for the constrained stress and strain fields. Noting that 



dXp 3 y >dx p 3 K ' r 

dlh-X 



dX 







hp (A19) 



we obtain 



du'i, d 



dX B dX 



4(l-i/) \r 



2 (1 - 2v) + (^\ | {2n a (h ■ X) - X a } 



+ 2{l-v) U 2 ) I 1 (r 2 



2 n-X 



-l)X a 



4(1- 



V) {-4 (1 - 2-) (£) - 4 (^) } { 2 "« • *) " *«} *f 



{2n a np - S a p} 




(A20) 



(l-2i/) + 



2n, 



( fl .x) 



I Xf3 

r [ r 



^) {2 (l-2i/)+(£)} {2h a hp-S a0 } 



1 - 2 



2 {n-X 



X n X 



r 2 J 1 \ r 2 



r r 



[^1 



2 n- X 



1>8. 



a/3 



13 



Thus the expression for t c a ^ = 2 f 



du " a dug 
dXp + dX a 



IS 



r 



r 



l^) (2(1-21/)+ (£) 



1 - 2 



2 n-X 



r 2 J I \ r 2 



I 2 ( fi ' £) x„x 



(A21) 



+4( I 



2 n-X 



1>« 



a/3 



allowing us to write the final expression for the constrained strain in the matrix 
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The trace of e c a g is not zero in the elastic medium and is given by 



'/'/ 



4(1-!/) 





2 m.\ -i 



1 - 2 



2 U-X 



2 n-X 



(A23) 



4(1-1/) 



(1 - 2i/) + 



a 2 \ fa 2 ^ " 2 



1-2^ (*\\ 2 Jt^L 1 



2 ^ -1 



(A24) 



We are now in a position to calculate the constrained stress in the elastic medium due to the deformed inclusion. 
It is given by the expression 



£v _ £e* 



• 4(l-i/ 2 ) 



£ 



vt a 



2\ 2 h-X 



(1 - v 2 ) \r 2 



- 1 ) 6 afj (A25) 



where 



is the expression inside the square brackets in Eq. IA22I 



2. Constrained displacement field - Cartesian components 



It proves useful to have the explicit cartesian components of the displacement field for computational and graphical 
purposes. If we consider the unit-vector h making an angle of <p with the positive direction of the x-axis, then the 



cartesian components of equation IA18I arc 



x 4(1-1/) V^ 2 



2(1-21/)+ 3 



[2 cos (x cos 4> + y sin < 



+ 



2(l-i/) 



+ 



x 4(l-i/) V^ 2 

'l 



2(1 - 2v) 



2 (x cos + y sin < 
,2 \ i 



[a; cos 2(f> + y sin 20] 



2(l-i/) 



x — y ) cos 20 + 2xy sin 20 



4(1-1/) \r 



+ 



2 (1 - 2v) + - 



[2 sin (x cos + y sin 0) — y] 



2(1-1/) Vr 



2 (x cos + y sin < 



a 4(1-1/) 

2 



2(1 - 2z/) + 



2(1-1/) V^ 2 



r 

x 2 — y 2 ) cos 20 + 2xy sin 20 



: sin 20 — y cos 20] 



Appendix B: Calculation of the energy of N quadrupoles 
Eq. (HHJ) can be re- written using e a p = 1/2 (it Q ,/3 + ws, a ) as 

i=l 17 v o ^ V-.L,i = i 

Using the symmetry of the stress tensor, we obtain 

JV 

2 i=i Jv$ 

We also have the identity 



i 



if there are no body forces. Thus we can write Eq. IB2l as 



dV 



Using Gauss's theorem to convert these volume integrals into area integrals, we obtain 

N N 

1 — 1 J SI J SI ^ JSrr, 
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where and n^°°^ are unit normal vectors both pointing outwards respectively from the inclusion volume Vq and 
the matrix boundary. Eq. IB5I can be rewritten as follows 



i r i N r 

s(oo) i=l J s o 



p _ - ^.(°°) c ( 00 ) 



2 «/3 c /3 7 
1 



AT 



4=1 JS, 



Thus we can write using the expressions earlier written 
down in Eq. f| 15[) 



iV 



i=l 
AT 

-S ) (^) = -i?+E^ i) ( J? ) ( B7 ) 

i=l 
N 

i m} (x)=u^(x) + Y,u^ (x) 



where [Xj indicates the constrained strain at lo- 
cation X in the matrix due to the eshelby labeled with 
the index i etc. We also have for locations X inside the 
inclusions 



(c,i) 



where e^-*^ is the eigenstrain of the ith Eshelby and so on. 
Note that in the expression for the strain in the inclusion 
given by Eq. (|B8[) we have removed the eigenstrain from 



the constrained strain e^g — e *g^ leaving only the elastic 
contribution in order to calculate correctly the clastic 
contribution to the energy. Using these expressions, the 
elastic energy of the system can be written from Eq. (|B6[) 



1 



^ — 2 a/3 € 0a V 
N 



(B9) 



i=l J13 o 



a/3 "73 w af3 "0 J "a 



Since the traction force has to be continuous at the in- 
clusion boundary (Newton's third law), we have 

^aB^a = ^aB^a at tne inclusion boundary (BIO) 



which gives us from Eq. (IB9 



iV 



i=l • /a o 



We also have from Eqs. (|B7[) and 



(i) (m) _ 

U /3 W /3 ~ e Bu 



(B12) 



On plugging this expression into Eq. (|B11[) gives us finally 



p _ 1 (oo) (oo) v 1 f (i) />(<),(*>*) v- jq 

A - 2 a/3 ^ a _ 2 ^ / w <*P a P" v 

i=l "'■So' 

^ (B13) 



N 



^ ^ — 9 ol B € Ba V 9 Zj ^0 € Ba a/3 



where a l aj3 = (1/Vq) J yi o-^dV. Using the expression for 
a aB f rom Eq. (|B8|) . we obtain 



W - + E - (CJ) (^) + <T - - ( *< l) (B14) 



Eq. (|E>14|) is a far field approximation that assumes that 
Rij 3> a. As Rij — > a clearly the spatial integrals con- 

tibuting to a^i must be computed explicitly and cannot 
be replaced by the single distance Rij between the centers 
of the Eshelby inclusions i and j. 

Using expression (|B13|1 we obtain 
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/ N \ N N 

A - 2 a? (ia ~ 2 Q/3 I ^ ^ I f3a afj 2 ^ /3q q/3 

\i=l / i=l i=l 



where all these terms are defined in Eqs. (|20| - ([23|) . 



(B15) 



Appendix C: The final form of Eh 



J 



We can also explicitly write Ei nc showing its linear dependence on the eigenstrain by using Eq. (|23[) . This gives us 



Eino = ~^*£ [{^hf ~ 6j) (ft;) + (2n«n« - tj) <r#W 



(Cl) 



Plugging Eq. (|A25|) into the above equation, we find that the term inside the square braces in Eq. (|C1[) can be written 



as: 



4(l-i/ 2 ) 



' ~i \rii\ /aWvW) i(')y('j)\ v('i! 
n • A ■> \ I n a JLp rip A Q \ A Q A^ 



ft 



ft 



(ft,) 2 



a 



+ 4 



ft 



ft 



2(l-2i/)+(^_J |{2n«n«-^} 
J | 1 



(ft,) 2 



a 

ft i 



ft 



ft 



ft 



cft^ ra^ 



ft 



1 - 



a 



21 2 (nW -J?W) 



(ft,) 2 



-u, 



aB 



(C2) 



where XW' indicates the vector joining the centers of the eshelby pair labeled as i and j, and (i O j) in Eq. (|C2[) 
represents the term obtained by exchanging i and j. 
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In order to simplify Eq. (|C2[) , we need the following identities 

(2r$nf - 5 af3 



flW.jfttA (f®x™ hfx^' 



8n w ■ h {3) 



Ri 



Ri 



Ri 



y(«i) y('j) 



2 



J? 




Ri 



{i «-> j) 
( fiW • (*•?') \ / n^') ■ 



Ri 



Ri 



■ h (j) 



Ri 



{R*j) 2 



n(«) . jf («') \ / fiO') . jf \ 



I?, 



(2h$nf - <5^) <S Q(9 + (i j) = 



(C3) 
(C4) 

(C5) 



(C6) 
(C7) 



Using these identities, we can write the final expression for the interaction energy in the form shown in Eq. 

I 
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